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Abstract 

Bilayer graphene has the very interesting property of an energy gap tunable with the vertical electric field. We 
propose an analytical model for a bilayer-graphene field-effect transistor, suitable for exploring the design parameter 
space in order to design a device structure with promising performance in terms of transistor operation. Our model, 
based on the effective mass approximation and ballistic transport assumptions, takes into account bilayer-graphene 
tunable gap and self-polarization, and includes all band-to-band tunneling current components, which are shown to 
represent the major limitation to transistor operation, because the achievable energy gap is not sufficient to obtain a 
large I n/I « ratio. 
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I. Introduction 

The progress of CMOS technology, with the pace foreseen by the International Technology Semiconductor 
Roadmap (ITRS) [1], cannot be based only on the capability to scale down device dimensions, but requires the 
introduction of new device architectures [2] and new materials for the channel, the gate stack and the contacts. This 
trend has already emerged for the recent technology nodes, and will hold - probably requiring more aggressive 
innovations - for devices at the end of the Roadmap. 

In the last decade carbon allotropes have attracted the attention of the scientific community, first with carbon 
nanotubes [3] and, since its isolation in 2004, with graphene [4], which has shown unique electronic [5] and 
physical properties [6], such as unconventional integer quantum Hall effect [7], [8], high carrier mobility [4] at 
room temperature, and potential for a wide range of applications [9], [10], [11], like nanoribbon FETs [12]. Despite 
graphene is a zero gap material, an energy gap can be engineered by "rolling" it in carbon nanotubes [13] or by the 
definition of lateral confinement like in graphene nanoribbons [14]. However, theoretical [15] and experimental [16] 
works have shown that significant gap in nanoribbons is obtained for widths close to 1-2 nm, which are prohibitive 
for fabrication technology on the scale of integrated circuits, at least in the medium term. 

Recently, theoretical models [17], [18], [19] and experiments [20] have shown that bilayer graphene has the 
interesting property of an energy gap tunable with an applied vertical electric field. Anyway, the largest attainable gap 
is of few hundreds of meV, which make its use questionable for nanoelectronics applications: limits and potentials 
of bilayer graphene still have to be shown. 

From this point of view, device simulations can greatly help in assessing device performance. Bilayer-graphene 
FETs (BG-FETs) have been compared against monolayer FETs, by means of the effective mass approximation [21] 
and Monte Carlo simulations [22] in the ballistic limit, showing really poor potential as compared to ITRS 
requirements [1]. These approaches, however, did not take into account some of the main specific and important 
properties of bilayer graphene, such as the possibility of tuning the band gap and the dispersion relation with the 
vertical electric field, and dielectric polarization in the direction perpendicular to the 2D sheet. Such problems have 
been overcome in Ref. [23], using a real space Tight-Binding approach. However, for the limited set of device 
structures considered, the small band gap does not allow a proper on and off switching of the transistor. 

One limitation of detailed physical simulations is that, despite their accuracy, they are typically too demanding 
from a computational point of view for a complete investigation of device potential. Analytical approaches could 
help in this case. One example has been proposed in Ref. [24], but it has serious drawbacks, because it completely 
neglects band-to-band tunneling and the dependence of the effective mass on the vertical electric field, providing a 
unrealistic optimistic picture of the achievable performance. 

In this work, we have developed a semi-analytical model for a bilayer-graphene FET with two gates to study 
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the possibility of realizing an FET by tuning the gap with a vertical electric field. The model has been validated 
through comparison with results obtained by means of a full 3D atomistic Poisson-Schrodinger solver, showing 
good agreement in the applied bias range [23], [25]. Interband tunneling proves to be the main limiting factor in 
device operation, as demonstrated by the device analysis performed in the parameter space. 

II. Model 

In this section we provide a detailed description of the developed model, which is based both on a top of the 
barrier model [26] and on the calculation of all the interband tunneling components. In particular we adopt the 
ballistic transport and the effective mass approximation, whose main electrical quantities, such as the effective mass 
and the energy gap, have been extracted from the energy bands obtained from a ^-orbital Tight Binding (TB) 
Hamiltonian. Since we want to address long channel devices, short channel effects have been completely neglected, 
as well as inelastic scattering mechanisms, which are expected to be negligible in this kind of material [6], With 
respect to more accurate atomistic models, the followed approach may underestimate the actual concentration of 
carriers in the channel, especially for large drain-to-source (Vds) and gate voltages (Vgs), when parabolic band 
misses to match the exact dispersion relation. We however believe that the developed model represents a good 
trade-off between accuracy and speed. 

A. Effective mass approximation 

In order to proceed with the definition of an analytical model based on the effective mass approximation, we 
first need an expression for the energy bands of bilayer graphene. The top view of the bilayer-graphene lattice 
structure with carbon-carbon distance a = 1.44 A is shown in Fig. \V[a): Al-Bl atoms lay on the top layer, while 
A2-B2 on the bottom layer. The energy dispersion relation can be computed by means of a p 2 -Tight Binding (TB) 
Hamiltonian [27] considering two layers of graphene coupled in correspondence of the overlaying atoms Al and 
A2. 

The energy dispersion relation reads [17] : 

£ (k ) = l h±Hi ± (1 ) 



l/(k)| 2 + ^ + t f ± l^(u 2 + ti)\fm 2 + t^ 

where U\ and U2 are the potential energies on the first and second layer, respectively, U = U\ — U2, ij_=-0.35 eV 
is the inter-layer hopping parameter [17], k = k x k x + k y k y and [27]: 



/(k) = te 



2 ( kyOVl] e _ i3kxa/2 



(2) 



which is the well known off-diagonal element of the 2x2 graphene p z -Hamiltonian, where t is the in-plane hopping 
parameter (£=-2.7 eV). In Fig.[Tfb) the band diagram for U — 0.5 eV is shown. As can be seen, bilayer graphene has 
four bands, symmetric with respect to the coordinate axis. For large U, the "mexican-hat" behavior in correspondence 
of the band minima can be observed, as detailed in Fig. QI C )- 

Let us now consider the third band (Fig. Q2b)), which corresponds to the conduction band (same considerations 
follow for the valence band, i.e. second band) and apply a parabolic band approximation in correspondence of the 
j, which reads [17]: 



mm y ip + 1\ 2v F n 

The dispersion relation can now be expressed as [17] 



E(k) = 5p + JL(|k| - <w„) 2 + (4) 



where 

,_^(^ 2 +ii) 3 /2 1 _ Ut± 



Eqap r—r\ n-i (5) 
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Fig. 1. a) Real space lattice structure of bilayer graphene. The bilayer consists of two coupled hexagonal lattices with inequivalent sites Al, 
Bl and A2, B2 in the first and in the second sheet, respectively, arranged according to Bernal (A2-B1) stacking, b) Tight-Binding band structure 
of bilayer graphene for U=Ui-U2=0.5 eV. c) Detail of the band structure in correspondence of band minimum k m i n : K is the Dirac point. 



vp — ^£ is the Fermi velocity and h is the reduced Planck's constant. 

As can be observed in (0, the effective mass m* has a singularity for U — 0, which is clearly unphysical. In 
order to avoid such an issue, energy bands in the range U £ [0, 0.14] have been fitted with the parabolic expression 
in (O, within an energy range of 2/cbT from the band minimum (where ks is the Boltzmann constant and T is the 
room temperature), and using m* as a fitting parameter. In Figs. 0a), Etb), we show, for two different inter-layer 
potential energies (U=Q eV and U=Q.l eV), the TB energy bands as well as the parabolic bands exploiting the 
analytical expression in (0 and the fitted values for m*, respectively. As can be seen, the fitted effective mass 
manages to better match the TB band in the specified energy range. In Fig. |2c), we show the fitted effective mass 
for different U . In particular, for U < 0.14 eV, m* can be expressed as: 

m* = 0.09E/ + 0.043, (6) 

while for larger values eq. (0 recovers. 



B. Electrostatics 

Once obtained the expression for m*, the electron concentration n can be expressed as: 



v 
2 



/•+0O 

/ D(E)[f(E-E FS )+f(E-E FD )]dE, (7) 
Je 
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Fig. 2. Comparison between energy dispersions obtained by means of the analytical effective mass (dashed-dotted line), fitted effective mass 
(solid line) and TB Hamiltonian (circle) for an inter-layer potential equal to a) (7=0 eV and b) [7=0.1 eV. c) Analytical and fitted relative 
effective mass as a function of the inter-layer potential U. m e is the free electron mass, d) ct co7l d{U) and ce va i(U) as a function of the 
inter-layer potential U. 



where / is the Fermi-Dirac occupation factor, Eps and E FB are the Fermi energies of the source and drain, 
respectively, and v=2 is band degeneracy. D(E) is the total density of states per unit area (for the complete 
calculation see the Appendix), which reads: 



D(E) = -L 
2nn 



2m* 



2m* 



E-E r 



where E c is the conduction band edge. If we define: 



fn(Ef) 



m 



k B T\n 



1 + exp 



E r - E 



k mm y/2m*k B T 
2^h 



/2 



k B T 

E c - Ef 
k B T 



where F 1 / 2 is the Fermi-Dirac integral of order 1/2, the electron concentration reads: 

n=[f n (E FS ) + f n (E FD )]. 
Analogous considerations can be made for the hole concentration p, which reads: 

p=[f P (E FS )+f p (E FD )}, 



where 



k min \/2m*kBT 
2ttH 



1 + cxp 

Fl/2 



Ef — E v 
k B T 

Ef ~ Ey 

k B T 



(8) 



(9) 



(10) 



(11) 



(12) 
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and E v is the valence band edge. 

Once n and p are computed, attention has to be posed on how charge distributes on the two layers i.e. on dielectric 
polarization. To this purpose, we have numerically extracted from TB simulations a va i(U) and a con d(U), that 
represent the fraction of the total states in the valence band and of electrons in the conduction band, respectively, 
on layer 1 [23]. We computed a CO nd(f ) f° r a particular bias (U% = — U 2 = U/2 and Ep = eV) and made 
the assumption that its dependence on the bias can be neglected. As far as a va i is concerned, we assumed in our 
considered bias range, that all electron states in the valence band are fully occupied and therefore f(E) = 1. 
Fig. |2jd), shows a con d(U) and a va i(U) as a function of the inter-layer potential U. The charge density pj per unit 
area on layer j (J =1,2) is expressed as the sum of the polarization charge, electrons and holes and finally reads: 

Pl (U) = q{[l-2a val (U)]N tot 

-n [1 - 2a cond (U)} + pa cond (U)} ; (13) 
p 2 (U) = q{[2a^{U)-l]N tot 
(U)+p[l- 

■^Q^cond 

where q is the electron charge and N tot is the concentration of ions per unit area. 

The considered device structure is a double-gate FET embedded in Si02- The bilayer graphene inter-layer distance 
d is equal to 0.35 nm, while two different oxide thicknesses t\ and t 2 have been considered (Fig. 0a)). An air 
interface between bilayer graphene and oxide has also been taken into account (t sp —0.5 nm) [28]. For such a system, 



ei to 



c 2 = 



*2 4. tl 
£2 Cf) 



we can define an equivalent capacitance circuit as in Fig. Ob), where Cq=%, C\= 
and e±= e 2 = 3.9eo, while eo = 8.85 x 10 -12 F/m. Vr g and Vsg are the top gate and back gate voltage respectively, 
V\ = = ^ L and V 2 = ~^ 2 -- In Fig. Oc), the flat band diagram along the transverse direction (y axis) is shown. 
Metal work functions for the back gate and top gate are equal to 4.1 eV [ < &s s = < l ) Tg=4. 1 eV], while the graphene 
work function ($ gra ) is equal to 4.5 eV [29]. E^Tg, Efs 9 are the Fermi level of the top and of the back gate, 
respectively. 

The conduction band edge inserted in eq. (O, can be expressed as: 

E c = + E FBg $ gra + + ^p. (14) 

Applying the Gauss theorem, we obtain the following expression: 

C^Vrg -V 1 ) + (V 2 - ViJCfa = -Pi 



Co(V 1 -V 2 ) + (V Bg -V 2 )C 2 = -p 2 . (15) 
Eqs. (fT3l l and ( fTBI l are then solved self-consistently till convergence on V\ and V 2 is achieved. 

C. Current 

Drain-to-source (Jtot) current is computed at the end of the self-consistent scheme. As depicted in Fig. 0d), 
Jtot consists of three different components: the first is due to the thermionic current J t h over the barrier [26], 
whereas the second (Jts) an d the third (Jtd) to band-to-band tunneling. In the same picture, we sketch the 
conduction band edge Ecs (Ecd) and the valence band edge Eys {Evd) at the source (drain). Assuming 
reflectionless contacts, the thermionic current is due to electrons injected from the source with positive velocity 
v x > and to electrons injected from the drain with v x < 0: 

f +00 



Jth — 



r+00 r /■ c^rp 

^L dk "\Loiz f{E - EFs)dk - 

/ ^rf(E-E FD )dk x 

Jk< ok x 



(16) 



where E = E c + (|k| — k min ) 2 , v x — is the group velocity and k x (k x ) is the wavevector range for 

which v x > (v x < 0). For the complete derivation see the Appendix. 
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Fig. 3. a) Sketch of the considered bilayer-graphene FET: d is the inter-layer distance, t\ and t2 are the top and back oxide thicknesses, b) 
Equivalent capacitance circuit of the simulated device, c) Flat band diagram of the BG-FET along the y direction, d) Conduction and valence 
band edge profiles in the longitudinal direction; we assume that deep in the source and drain regions, the electric field induced by the gate 
vanishes and the gap gradually reduces to zero. 



Let us now discuss the band-to-band tunneling current due to the barrier at source(drain) contact, which reads: 

Jn = 2 [ I q J^l^- Ti (k y )[f(E-E FS ) 
Jk v Jk> 2tH h dk x 

-f(E- E FD ))dk x dk y i = S,D, (17) 

where S refers to the source and D to the drain, while Ti(k y ) is the transmission coefficient at the different reservoirs. 
The key issue in computing ( fTTI i is the definition of an expression for Ti(k y ), which accounts for band-to-band 
tunneling process. 

We have assumed a non charge-neutrality region of fixed width Ax at the contact/channel interface and an electric 
field £i=(E c — Epi) / (qAx) with i = S,D. For what concern the Jts term, electrons emitted with electrochemical 
potential Eps see two triangular barriers, one at the source junction and one in correspondence of the drain 
(Fig. [3d)), whose heights are equal to E gap and width Wi = E gap /(q£i). Assuming the same Ax for both source 
and drain junctions, the drain barrier is transparent with respect to the source barrier, since, for large Vds, the 
electric field at the source is smaller than the electric field at the drain barrier: Ts(k y ) is therefore essentially given 
by the source junction barrier. Same considerations follow for the other band-to-band tunneling current component 
Jtd, flowing only through the drain-channel contact. In this case Sd = E c — Efd = Ec ~ E ^ s x +qVlDS . 
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Assuming the WKB approximation, the transmission coefficient can be expressed as: 

Ti{k y ) = e~ 2 S ^ Im{kxWx , i = S,D, (18) 
where Im{k x } is the imaginary part of k x and is obtained from: 

h 2 

(|k| - k min ) = q£iX - E gap , i = S,D. (19) 
Finally, Jti is computed performing the integral ( fTTT i numerically. 

III. Exploration of the design space 

In order to validate our model, we have first compared analytical results with those obtained by means of 
numerical NEGF Tight Binding simulations [25], considering a test structure with t\ =ti =1.5 nm, t sp = 0.5 nm, 
Qgra = §g = $s fl =4.1 eV, Vds=0.1 V and VB g =0 V. In Fig. BJa)-(b) the electron concentrations (pi, pi) and 




Fig. 4. Comparison between analytical and numerical simulation of a) pi and p2 and b) Vi and V2 as a function of Vr g - Vb 9 =0 V and 
Vds=0 V. c) Energy gap as a function of top gate voltage, with Vbq = V. 



the electrostatic potentials (Vi, V%) on layer 1 and 2 are shown, as a function of Vr g , for Vds=0 V and Vb 9 =0 V. 
As can be seen, results are in good agreement. Some discrepancies however occur for larger Vds (Vds > 0.2 V), 
where the parabolic band approximation misses to reproduce band behavior for large k y . In Fig. |4|c) the energy 
gap is plotted as a function of Vr g - As can be seen, even for large Vr g , the biggest attainable E gap is close to 
0.15 eV. 

Let us now consider the different contributions of the three current components (Jth, Jts an d Jtd) to the total 
current Jtot (Fig- Hd)). For each of these components, we can define a sort of threshold voltage, above which 
their contribution is not negligible. In particular, J t h starts to be relevant as soon as E c ~ Eps- We then define 
Vth as the Vr g for which E c — Eps- Similarly, interband current Jts is not zero when E v >Ecs so we define 
Vts the top-gate voltage for which E v —Ecs- Finally, Jtd is not zero in the energy range Ecd < E < Evs'- we 
define Vtd anc ' ^rb tne top-gate voltages for which Evs=Ecd', thanks to these definitions, we can qualitatively 
evaluate current contribution by observing the band structure. 

Our goal is indeed to obtain the largest value for the I on /I ff ratio, and this is only possible if the band-to-band 
component of the current is suppressed. We have considered three different solutions to accomplish this task: by 
varying the back gate oxide (t 2 ), by varying the Eps ~ Ecs ° r Epp, — Ecd difference, or by simply varying the 
back gate voltage. If otherwise specified, Ax = 0.7 nm, as obtained from TB simulations of an abrupt junction 
with the same doping of the considered BG-FET. In Fig. |3a)-(b)-(c) the above-defined thresholds are shown for 
the three considered cases. 

As shown in Fig. |2a), back gate oxide thickness has no effect in our case, since the top layer screens the electric 
field induced by the top gate, as can also be seen from Fig. H|b), where V2 remains almost constant. Fig. |5}b) 
shows thresholds as a function of (Eps — Ecs), an d therefore as a function of dopant concentration. We observe 
that for Efs — Ecs — 1 eV, Vtd = ^rb> so triat Jtd is practically eliminated, while Jts increases since Vts 



September 25, 2009 



DRAFT 



8 




Fig. 5. a) Thresholds as a function of the oxide thickness for Vp>s = 0.5 V, Vbg = V, Eps — Ecg=0.9 eV. b) Thresholds as a function 
of Eps — Ecs f° r t2=1.5 nm, Vds= 0.5 V, Vbg = ® V. c) Thresholds as a function of Vbq for t2=l .5 nm, Vos— 0.5 V and Eps — Eqs- 
1 eV. d) Total current Jtot and its three components (Jts, Jtd< Jth) f° r ^DS-0-1 V, Vs9=0 V and Epg — Eqs= 0.5 eV . Thresholds 
are shown along the coordinate axis, e) Total current for £a=1.5 nm, Vp>s- 0.5 V, Epg — Eq§= 1 eV and Vgc=0 V and different Ax. 



becomes larger. In Fig.|5jc), we show Vth and Vrs, for Eps — Ecs = 1 eV^, as a function of Vb^- Unfortunately 
the two curves have the same behavior, so that Vrs cannot be reduced to values smaller than Vth, or -in other 
words- we cannot suppress current due to interband tunneling at source contact. 

We have then computed the transfer characteristics for Vbq = 01^, t\ = t 2 =1.5 nm, Eps—Ecs=^ eV. In Fig. [5J e ) 
Jtot is shown. As can be seen, poor I on /I s ratio can be obtained since band-to-band tunneling at source contact 
is too high as also observed in graphene FET [30]. Reducing £, i.e. T(k y ), could lead to a reduction of Jts an d 
consequently to an improvement of the I on /I s ratio. As can be seen in Fig.[5]e), an improved I 0B ./I s is obtained 
increasing Aa; to 5-10 nm, but it is still lower than the ITRS requirements (10 4 ) for digital circuits. In Fig. [6] we 
also sketch the simulated band edges for three different cases: a) when tunneling is negligible ( Jtot — Jt% for 
Vtq = 2 V), b) when tunneling weakly affect the total current (Jtot — 10(Jts + Jtd) for Vr g — 0.4 V) and c) 
when tunneling represents the predominant component (Jtot — (Jts + Jtd) for Vr g = —2 V). 

IV. Conclusion 

We have developed an analytical model for bilayer-graphene field effect transistors, suitable for the exploration 
of the design parameter space. The model is based on some simplifying assumptions, such as the effective mass 
approximation, but includes all the relevant physics of bilayer graphene. First and foremost, it includes the tunable 
gap of bilayer graphene with the vertical electric field, which is exploited in order to induce the largest gap, when 
the device is in the off state. It also fully includes polarization of bilayer graphene in response to a vertical electric 
field. As far as transport is concerned, it includes the thermionic current components and all interband tunneling 



September 25, 2009 



DRAFT 



9 




Fig. 6. Simulated band edges when: a) Tunnel current weakly affects Jtot\ t>) Jts + Jtd is one order of magnitude smaller than the total 
current; c) Inter band tunneling current is the dominant component. 



components, which are the main limiting factor in achieving a large I on /I g ratio. Significant aspects of the model 
have been validated through comparisons with numerical TB NEGF simulations. 

Due to the small computational requirements, we have been able to explore the parameter design space of bilayer- 
graphene FETs in order to maximize the I on /Iotf ratio. Despite applied vertical field manages to induce an energy 
gap of the order of one hundred meV, band-to-band tunneling greatly affects device performance, limiting its use 
for device applications. A larger gap must be induced to make bilayer graphene a useful channel material for digital 
applications, probably by combining different options, such as using bilayer graphene in addition to limited lateral 
confinement, stress, or doping. 



Appendix I 

A. Density of states 

The total density of states can be computed as follow, performing the integral over the first Brillouin zone (BZ) 

2 r „/„ h 2 



D ( E ) = 7^2 / 5 [ E o + ^r(l k l - k ^f - E ) 2nkdk. (20) 
Jbz V 2m 



If we apply the following property of the delta function 



\f'(Xn)\ ' 

where x n are the zeroes of the function f(x), eq. (1201 ) reads: 



D(E) = l^f J f S(r-VE^) 
V ' it h 2^E - E c J BZ , V v V 



^fr + kmm ) d r 



2-Kh 

where r=-J=(\k\-k min ). 



1 (2m* I 2m* 

" ruin I ; \ z ' z '/ 
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B. Thermionic Current 

In order to derive the expression of the thermionic current, we have first to compute the group velocity v x , which 
reads: 

1 dE(k x , k y ) -. a / F *2 

h dk x 



1 d fE gap h f | k ,_. ,2 



h dk x 

m* 



2 2m 

k 

>' in, 



1 



(23) 



Replacing ( |23l in ( fT6] l we obtain: 



Jth — 



—qh 
TT 2 m* 



dk*. 



k x [l- 



f(E - E FS )dk x 



k x 1 



f(E-E FD )dk x 



(24) 



In order to remove the singularity in eq. d24i i for k x = k y = 0, we can use cylindrical coordinates, i.e. k x =kcos8, 
k y =ksin9 and |k|=fc. In this representation, the condition v x > translates in: 

kcosO < 



kcosO > 
k k m i n ^> 



k kmiji, <C 



The integral (1241) becomes: 



Jth 



2qh 



k(k- k min ) [f s {E, k) - f D (E, k)] dk. 



(25) 



C. Transmission coefficient 

The tunneling transmission probability T(k y ) has been computed through the WKB approximation. The |/m{fc x }| 
in eq. (TTsT t is computed from the energy dispersion relation as follows: from (T% we can write 



h 2 — h 



2,771 . Eq a p . 

= -jp-(^-— )• 



Defining 



and inserting eq. d2~7l > in eq. 



2m ; 



ii 

we obtain: 

1 kmin — 



> 



2 £^ 



2i(i(x)k ri 



(26) 

(27) 

(28) 
(29) 
(30) 

(31) 

withC = -k 2 mm + p 2 (x)+k 2 y . 
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which reads: 

k x = (^min ~ Pi X ) 

If we expressed k x as k x =a + ib {a, b 6 5ft), k 2 reads: 



kl = a 2 - b 2 + 2iab. 



By comparing eq. ( 1291 and eq. (13 PI . |/m{A: K }| simply reads: 

\Im{k x }\ = 



c + Vc 2 + W^Lr 
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